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Oj' Abstract 

The complexity of the atmosphere endows it with the property of tur- 
bulence by virtue of which, wind speed variations in the atmospheric 
^ ■ boundary layer (ABL) exhibit highly irregular fluctuations that persist 

OA ' over a wide range of temporal and spatial scales. Despite the large and sig- 

Cn , niflcant body of work on microscale turbulence, understanding the statis- 

^P ■ tics of atmospheric wind speed variations has proved to be elusive and 

challenging. Knowledge about the nature of wind speed at ABL has far 
reaching impact on several fields of research such as meteorology, hy- 
drology, agriculture, pollutant dispersion, and more importantly wind en- 
ergy generation. In the present study, temporal wind speed records from 



Y^ . twenty eight stations distributed through out the state of North Dakota 

, ,_( ' (ND, USA), (~ 70,000 square-miles) and spanning a period of nearly eight 

^ , years are analyzed. We show that these records exhibit a characteristic 



broad multifractal spectrum irrespective of the geographical location and 
topography. The rapid progression of air masses with distinct qualitative 
characteristics originating from Polar regions. Gulf of Mexico and North- 
ern Pacific account for irregular changes in the local weather system in 
ND. We hypothesize that one of the primary reasons for the observed 
multifractal structure could be the irregular recurrence and confluence of 



j^ ' these three air masses. 



1 Introduction 



Atmospheric phenomena are accompanied by variations at spatial and temporal 
scales. In the present study, qualitative aspects of temporal wind speed data 
recorded at an altitude of 10 ft from the earth's surface are discussed. Such 



recordings fall under the ABL, which is the region 1-2 km from the earths sur- 
face ^. Flows in the ABL, which are generally known to be turbulent, are 
influenced by a number of factors including shearing stresses, convective insta- 
bilities, surface friction and topography PH]. The study of laboratory scale 
turbulent velocity fields has received a lot of attention in the past (see [3] for a 
summary). A. N. Kolmogorov 01 13], (K41) proposed a similarity theory where 
energy in the inertial sub-range is cascaded from the larger to smaller eddies 
under the assumption of local isotropy. For the same, K41 statistics is also 
termed as small-scale turbulence statistics. The seminal work of Kolmogorov 
encouraged researchers to investigate scaling properties of turbulent velocity 
fields using the concepts of fractals . Subsequent works of Parisi and Frisch 
|7], Meneveau and Srinivasan, [SJ |^ ^] provided a multifractal description of 
turbulent velocity fields. While there has been a precedence of scaling behavior 
in turbulence at microscopic scales ^ El [3 El El EH E] it is not necessary that 
such a scaling manifest itself at macroscopic scales, although there have been 
indications of "unified scaling" models of atmospheric dynamics, ^21- Several 
factors can significantly affect the behavior of a complex system such as ABL 
m d. Thus, an extension of these earlier findings |lElEl|3|llEl[Tni[n] 
to the present study is neither immediate, nor obvious. Attempts have also 
been made to simulate the behavior of the ABL |13[ I14|. However, there are 
implicit assumptions made in these studies and often, there can be significant 
discrepancies between simulated flows and the actual phenomenon when these 
assumptions are violated |3]. On the other hand, knowledge about the nature 
of wind speed at ABL has far reaching impact on several fields of research. In 
particular, the need to obtain accurate statistical descriptions of flows in the 
ABL from actual site recordings is both urgent and important, given its utility 
in the planning, design and efficient operation of wind turbines, ^Hl- There- 
fore, analysis of wind speed records based on numerical techniques is gaining 
importance in the recent years. In |16j . long term daily records of wind speed 
and direction were represented as a two dimensional random walk and the re- 
sults reinforce the important role that memory effects have on the dynamics 
of complex systems. In |17| . the short term memory of recorded surface wind 
speed records is utilized to build rn'th order Markov chain models, from which 
probabilistic forecasts of short time wind gusts are made. In ^j, the authors 
study the correlations in wind speed data sets over a span of 24 hours, using 
detrended fluctuation analysis (DFA), jEj and its extension, multifractal-DFA 
(MF-DFA) |2n] . Their studies show that the records display long range corre- 
lations with a fluctuation exponent of a ~ 1.1 along with a broad multifractal 
spectrum. In addition, they also suggest the need for detailed analysis of data 
sets from several other stations to ascertain whether such features are charac- 
teristic of wind speed records. In |23) it is shown that rare events such as wind 
gusts in wind speed data sets that are long range correlated are themselves long 
range correlated. In 22 , it is shown that surface layer wind speed records can 
be characterized by multiplicative cascade models with different scaling rela- 
tions in the microscale inertial range and the mesoscale. Our previous studies, 
pnj suggest that at short time scales, hourly average wind speed records are 



characterized by a scaling exponent a ^ 1.4 and at large time scales, by an ex- 
ponent of a ~ 0.7. A deeper examination of the data sets in J5B] using MF-DFA 
indicated that the records also admitted a broad multifractal spectrum under 
the assumption of a binomial multiplicative cascade model. Interestingly, scal- 
ing phenomena have also been found in fluctuations of meteorological variables 
that influence wind speed such as air humidity, |24) . temperature records and 
precipitation [2S|- In EHj it is observed that while temperature records dis- 
play long range correlations (a ~ 0.7), they do not display a broad multifractal 
spectrum. On the other hand, precipitation records display a very weak degree 
of correlation with (a ~ 0.5, [2S]). While it is difficult to directly relate the 
scaling results of these variables to that of wind speed, greater insight can be 
gained by analyzing data sets that are recorded over long spans from different 
meteorological stations. Motivated by these findings, we chose to investigate 
the temporal aspects of wind speed records dispersed over a wide geographical 
area. In the present study, we follow a systematic approach in determining the 
nature of the scaling of wind speed records recorded at an altitude of 10 ft across 
28 spatially separated locations spanning an area of approximately 70,000 sq.mi 
and recorded over a period of nearly 8 years in the state of North Dakota. 
As noted earlier, convective instabilities and topography can have a prominent 
impact of the flow in ABL. The air motion over North Dakota is governed by 
the flow of three distinct air masses with distinct qualities, namely from : (i) the 
polar regions to the north (cold and dry) (ii) the Gulf of Mexico to the south 
(warm and moist) and (iii) the northern pacific region (mild and dry). The 
rapid progression and interaction of these air masses results in the region being 
subject to considerable climactic variability. These in turn can have a significant 
impact on the convective instabilities which governs the flow in ABL. On the 
other hand, the topography of the region has sharp contrasts on the eastern 
and western parts of the state because of their approximate separation by the 
boundary of continental glaciation. The eastern regions have a soft topography 
compared to the western region which comprises mostly of rugged bedrock. 
In the present study, we show that the qualitative characteristics of the wind 
speed records do not change across the spatially separated locations despite the 
contrasting topography. This leads us to hypothesize that the confluence of the 
air masses as opposed to the topography plays a primary role in governing the 
wind dynamics over ND. 

2 Methods 

Spectral analysis of stationary processes is related to correlations in it by the 
Wiener-Khinchin theorem, [2I| and has been used for detecting possible long- 
range correlations. In the present study, we observed broad-band power-law 
decay superimposed with peaks. This spectral signature was consistent across 
all the 28 stations (Fig. db)). Such power-law processes lack well-defined 
scales and have been attributed to self-organized criticality, intermittency, self- 
similarity |28l I2H| and multiscale randomness I^OI . Superimposed on the power- 



law spectrum, were two high frequency peaks which occur at ii ~ 24 and f2 = 12 
hours corresponding to diurnal and semi-diurnal cycles respectively. 
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Figure 1: (a) Temporal trace of hourly average wind speed record (miles/hour) 
at one of the representative stations (Baker IN, refer to Table 1 for details) 
over a period of nearly 8 years, (b) The corresponding power spectrum exhibits 
a power law decay of the form {S{f) ~ 1/f^). Superimposed on the power 
spectrum are prominent peaks which correspond to multiple sinusoidal trends. 

(c) Log-Log plots of the fluctuation function versus time scale, Fq{s) vs s for 
the moments q = -10(*), -6 (triangle ), -0.2(x) , 2 (.), 6(o) and q = lO(-h) 

(d) Multifractal spectrum of the record determined under the assumption of a 
binomial multiplicative cascade model. 



Power-law decay of the power-spectrum (Fig. ^b)) can provide cues to possi- 
ble long-range correlations, but, however, it is susceptible to trends and non- 
stationarities which are ubiquitous in recordings of natural phenomena. While 
several estimators [311 \'6'2\ have been proposed in the past for determining the 
scaling exponents from the given data, detrended fluctuation analysis (DFA), 
[T^ and its extension, generalized multifractal-DFA (MF-DFA),20, have been 
widely used to determine the nature of the scaling in data obtained from di- 
verse settings [201 ESI ISI EH ESI • In DFA, the scaling exponent for the given 
monofractal data is determined from least-squares fit to the log-log plot of the 
second-order fluctuation functions versus the time-scale, i.e. Fq{s) vs s where 
q = 2. For MF-DFA, the variation of the fluctuation function with time scale is 



determined for varying q {q^ 0). The superiority of DFA and MF-DFA to other 
estimators along with a complete description is discussed elsewhere |32j- DFA 
and MF-DFA procedures follow a differencing approach that can be useful in 
eliminating local trends JH|- However, recent studies have indicated the suscep- 
tibility of DFA and MF-DFA to higher order polynomial trends. Subsequently, 
DFA-n \67\ was proposed to eliminate polynomial trends up to order n — 1. In 
the present, study we have used polynomial detrending of order four. However, 
such an approach might not be sufficient to remove sinusoidal trends which can 
be periodic |3H| or quasiperiodic (see discussion in Appendix A). 

Data sets spanning a number of years, as discussed here, are susceptible to sea- 
sonal trends that can be periodic or quasiperiodic in nature. Such trends man- 
ifest themselves as peaks in the power spectrum and their power is a fraction 
of the broad-band noise. These trends can also introduce spurious crossovers 
as reflected by log-log plot oi Fq{s) vs s and prevent reliable estimation of the 
scaling exponent. Such crossovers indicate spurious existence of multiple scaling 
exponents at different scales and shift towards higher time scales with increasing 
order of polynomial detrending (27] . Thus, it is important to discern correlations 
that are an outcome of trends from those of the power-law noise. In a recent 
study, [3H|, a singular- value decomposition (SVD) based approach was proposed 
to minimize the effect of the various types of trends superimposed on long-range 
correlated noise. However, SVD is a linear transform and may be susceptible 
when linearity assumptions are violated. Therefore, we provide a qualitative 
argument to identify multifractality in wind speed records superimposed with 
periodic and/or quasiperiodic trends. Multifractality is reflected by a change in 
the slope of the log- log fluctuation plots with varying q with {q ^ 0) ^^. For 
a fixed q, one observes spurious crossovers in monofractals well as multifractal 
data sets superimposed sinusoidal trends. Thus nonlinearity or a crossover of 
the log-log plot for a fixed q might be due to trends as opposed to the existence of 
multiple scaling exponents at different scales. However, we show (see discussion 
under Appendix A) that the nature of log- log plot of Fq{s) vs s does not change 
with varying q for monofractal data superimposed with sinusoidal trends. How- 
ever, a marked change in the nature if the log- log plots Fq{s) vs s with varying 
q is observed for multifractal data superimposed with trends. Therefore, in the 
present study, the log-log plot of Fq{s) vs s with varying q is used as a qualita- 
tive description of multifractal structure in the given data. For the wind-speed 
records across the 28 stations, we found the peaks in the power spectrum to 
be consistent across all the 28 stations. Thus any effects due to the trend on 
the multifractal structure, we believe would be consistent across the 28 stations. 



3 Results 

MF-DFA was applied to the data sets recorded at the 28 stations. The log-log 
plots of the fluctuation functions (Fq^s) vs s) with varying moments q = -10, 



-6, -0.2, 2, 6, 10) using fourth order polynomial detrending for one of the rep- 
resentative records is shown in Fig^c). From Fig^Jc), it can be observed that 
the data sets exhibit different scaling behavior with varying q, characteristic of 
a multifractal process. This has to be contrasted with monofractal data whose 
scaling behavior is indifferent to the choice of q in the presence or absence of 
sinusoidal trends. 

To compute the q dependence of the scaling exponent h{q), we select the time 
scale in the range [2.2 - 3.7] where the scaling was more or less constant for 
a given q. Note that this corresponds to variations in a time span of [10^'^ — 
10^'^] ~ [158 — 5012] hours, for a given q. In this range, the slope of the 
fluctuation curves h{q) was calculated for every q and for ever station. The 
mean of the generalized exponents h{q) over all the twenty eight stations along 
with the standard deviation bars are shown in Fig. El^a). 
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Figure 2: (a) The mean (circle) and the standard deviation (vertical lines) of 
the generalized Hurst exponent, h{q) vs q across the 28 stations (b) The mean 
(circle) and the standard deviation (vertical lines) of the multifractal spectrum 
/(a) across the 28 stations, (c) Histogram of the multifractal widths (Aa) 
across the 28 stations (d) Histogram of the Hurst exponent h{2). 

It can be noted from Fig. I3a) that the slopes h{q) decrease as the moment q 
varies from negative to positive values which signifies that wind speed fluctua- 
tions are heterogeneous and thus, a range of exponents is necessary to completely 
describe their scaling properties. To capture this notion of multifractality, we 
estimate the classical Renyi exponents T{q) and the singularity spectrum jHOI 



under the assumption of binomial multiplicative process 0] 001 ED] (see Ap- 
pendix A for details). The singularity spectrum of one of the representative 
stations (Baker, IN) is shown in Fig^d) and its variation across the 28 sta- 
tions is shown in Fig. |2Ib). The fitting parameters a, b for the cascade model, 
the Hurst exponent h{2) and the multifractal width Aa for all the stations are 
summarized in Table 1. 

These results indicate multifractal scaling consistent across the stations. Earlier 
studies 02 Oni HOI have suggested the choice of random shuffled surrogates to 
rule out the possibility that the observed fractal structure is due to correla- 
tions as opposed to broad probability distribution function. The wind speeds 
in the present study follow a two-parameter asymmetric WeibuU distribution 
whose parameters were also similar across the 28 stations. MF-DFA on the 
random shuffle surrogates of the original records FigJ^d) indicate a scaling of 
the form Fq(s) ~ s°'^ with varying q, characteristic of random noise and loss 
of multifractal structure. The width of the multifractal spectrum was used to 
characterize the strength of multifractality. The histogram of the multifractal 
widths obtained across the 28 stations. Fig. |2Ic), was narrow with mean and 
standard deviation A(a) = (0.4866 ± 0.0599). The multifractal widths and the 
Hurst exponent h{2) across the twenty eight stations is also shown in Fig|31 
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Figure 3: The multifractal width for each of the 28 stations is indicated by 
circles. The Hurst exponent h(2) is indicated by upright triangles. The x-y 
plane represents the x and y coordinate in the of North Dakota. 



In the present study a systematic approach was used to determine possible scal- 
ing in the temporal wind-speed records over 28 spatially separated stations in 
the state of North Dakota. Despite the spatial expanse and contrasting topog- 
raphy the multifractal qualitative characteristics of the wind speed records as 



Table 1: Names and locations of the 28 recording stations. The fitting param- 
eters (a, b) of the cascade model, the Hurst exponent h{2) and the multifractal 
width (Aa) are also indicted. 



Station Number 


Name 


a 


b 


h{2) 


Aa 


1 


Baker IN 


0.513 


0.710 


0.702 


0.469 


2 


Beach 9S 


0.523 


0.76 


0.643 


0.539 


3 


Bottineau 14W 


0.525 


0.720 


0.699 


0.456 


4 


Carrington 4N 


0.521 


0.721 


0.638 


0.468 


5 


Dazcy 2E 


0.551 


0.722 


0.609 


0.390 


6 


Dickinson INW 


0.535 


0.715 


0.670 


0.418 


7 


Edgeley 4SW 


0.531 


0.752 


0.573 


0.510 


8 


Fargo INW 


0.505 


0.719 


0.641 


0.510 


9 


Forest River 7WSW 


0.556 


0.739 


0.670 


0.411 


10 


Grand Forks 3S 


0.527 


0.757 


0.646 


0.522 


11 


Hazcn 2W 


0.544 


0.733 


0.672 


0.430 


12 


Hettinger INW 


0.523 


0.748 


0.610 


0.516 


13 


HiUsboro 7SE 


0.522 


0.751 


0.620 


0.526 


14 


Jamestown 10 W 


0.547 


0.739 


0.636 


0.434 


15 


Langdon IE 


0.545 


0.710 


0.644 


0.381 


16 


Linton 5N 


0.531 


0.711 


0.623 


0.422 


17 


Minot 4S 


0.550 


0.706 


0.629 


0.459 


18 


Oakes 4S 


0.507 


0.755 


0.657 


0.574 


19 


Prosper 5NW 


0.504 


0.720 


0.671 


0.513 


20 


Mohall IW 


0.513 


0.753 


0.654 


0.554 


21 


Streeter 6NW 


0.521 


0.742 


0.675 


0.509 


22 


Turtle Lake 4N 


0.546 


0.729 


0.693 


0.418 


23 


Watford City IW 


0.524 


0.768 


0.631 


0.551 


24 


St. Thomas 2WSW 


0.550 


0.744 


0.633 


0.434 


25 


Sidney INW 


0.518 


0.772 


0.635 


0.593 


26 


Cavalier 5W 


0.506 


0.736 


0.651 


0.539 


27 


Wilhston 5SW 


0.514 


0.759 


0.635 


0.562 


28 


Robinson 3NNW 


0.514 


0.739 


0.620 


0.525 



reflected by singularity spectrum, were found to be similar. Thus multifractal- 
ity may be an invariant feature in describing the dynamics long-term motion 
of wind speed records in ABL over North Dakota. We also believe that the 
irregular recurrence and confluence of the air masses from Polar, Gulf of Mexico 
and the northern Pacific may play an important role in explaining the observed 
multifractal structure. 
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A Data Acquisition, MF-DFA algorithm and Dis- 
cussion 

A.l Data Acquisition 

The wind speed records at the 28 stations spanning nearly 8 years were obtained 
from part of the chmatological archives of the state of North Dakota. Stations 
were selected to represent the general climate of the surrounding area. Wind 
speeds were recorded by means of conventional cup type anemometers located 
at a height of 10 ft. The anemometers have a range of to 100 mph with an 
accuracy of ±0.25 mph. Wind speeds acquired every five seconds are averaged 
over a 10 minute interval to compute the 10 minute average wind speed. The 
10 minute average wind speeds are further averaged over a period of one hour 
to obtain the hourly average wind speed. 



A. 2 Multifractal Detrended Fluctuation Analysis (MF- 
DFA): 

MF-DFA, [201 a generalization of DFA has been shown to reliably extract more 
than one scaling exponent from a time series. A brief description of the algo- 
rithm is provided here for completeness. A detailed explanation can be found 
elsewhere [^Ol- Consider a time series {xk},k = 1...N. The MF-DFA algo- 
rithm consists of the following steps. 

1. The series {xk} is integrated to form the integrated series {yk} given by 

i—k 

y{k)^Y.^x{i)-x] k^l,...N (1) 

where x represents the average value. 

2. The series {yk} is divided in to rig non-overlapping boxes of equal lengths 
where Us = int{N/s). To accommodate the fact that some of the data 
points may be left out, the procedure is repeated from the other end of 
the data set |20j . 

3. The local polynomial trend yv with order v is fit to the data in each box, 
the corresponding variance is given by 

^'(^' ''^ " I ; E{y[^ - (" - "-)-'^ + ^] - y-(*) } (2) 

for V — 1, . . .Us. Polynomial detrending of order m is capable of eliminat- 
ing trends up to order m-1. |20 | 
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4. The qth order fluctuation function is calculated from averaging over all 
segments. 



^5:[F^(.,.r/^| (3) 

In general, the index q can take any real value except zero. 

5. Step 3 is repeated over various time scales s. The scaling of the fluctuation 
functions Fq{s) versus the time scale s is revealed by the log-log plot. 

6. The scaling behavior of the fluctuation functions is determined by analyz- 
ing the log-log plots Fq{s) versus s for each q. If the original series {xk} 
is power-law correlated, the fluctuation function will vary as 

Fqis) ~ s''('^) (4) 

The MF-DFA algorithm [201 was used to compute the multifractal fluctuation 
functions. The slopes of the fluctuation functions {h{q)) for each 
q = (—10,— 6,— 0.2,2,6,10) was estimated by linear regression over the time 
scale range s — [2.2,3.7]. The generalized Hurst exponents {h{q)) are related 
to the Renyi exponents T{q) by qh{q) = T{q) + 1. The multifractal spectrum 
f{ah) defined by, |40i a/i = d^^ f{<^h) = qoth - T{q). Under the assump- 
tion of a binomial multiplicative cascade model 00] the generalized exponents 

h{q) can be determined from h{q) = ^ in2 — T'^s parameters a and b for 

each station was determined using a nonlinear least squares fit of the preceding 
formula with those calculated numerically. Finally, the multifractal width was 
calculated using Aa — Aa^ = /i(— oo) — h{oo) = " in2^ ; EOl- 



A. 3 Discussion 

The power spectrum of the wind speed records considered in the present study 
exhibited a power law decay of the form {S{f) ^ 1//'^) superimposed with 
prominent peaks which corresponds to multiple sinusoidal trends. Such a be- 
havior is to expected on data sets spanning several years. The nature of the 
power spectral signature was consistent across all stations. This enables us to 
compare the nature of scaling across the 28 stations. Recent studies had indi- 
cated the susceptibility of MF-DFA to sinusoidal trends in the given data, [23 ■ 
Sinusoidal trends can give rise to spurious crossovers and nonlinearity in the 
log-log plot that indicate the existence of more than one scaling exponent. In 
a recent study |3H1, we had successfully used the singular- value decomposition 
to minimize the effect of offset, power-law, periodic and quasi-periodic trends. 
However, SVD is a linear transform and may be susceptible when linearity as- 
sumptions are violated. Estimating the SVD for large-embedding matrices is 
computationally challenging. Therefore, in the present study we opted for a 
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qualitative description of multifractal structure by inspecting the nature of the 
log-log plots of the fluctuation function versus time scale Fg{s) vs s with varying 
moments q. We show that the nature of the log-log plot does not show appre- 
ciable change with varying moments q for monofractal data superimposed with 
sinusoidal trends. However, a marked change in the nature of the log-log plot 
is observed for multifractal data superimposed with sinusoidal trends. More- 
over, the nature of the trends as reflected by the power spectrum is consistent 
across the 28 stations. This enables us to compare the multifractal description 
obtained across the stations. The effectiveness of the qualitative description is 
demonstrated with synthetic monofractal and multifractal data sets superim- 
posed with sinusoidal trends. 

A. 3.1 MF-DFA results of monofractal and multifractal data super- 
imposed with sinusoidal trends 

Consider a signal yi{n) consisting of monofractal data si{n) superimposed with 
a sinusoidal trend ti{n). Let 2/2(^1) be a signal consisting of multifractal data 
S2{n) superimposed with sinusoidal trend t2{n). The trends are described by, 

ti{n) = Ai sin(27rn/Ti) 4- A2 sin(27rn/T2) + A3 sin(27rn/T3), n = 1 . . . TVi 
t2{n) = Bi sin(27rn/rib) + B2 sm{2Tm/T2b),n = \ . . . N2 

The signals are given by, yi{n) — Si{n)+ti{n), i — 1,2 where si(ri) is monofractal 
data with a = 0.9 and S2(n) is multifractal data is that of internet log traffic, 
|44| . The dominant spectral peaks Fig^a) and Fig^b) reflect the presence 
of these trends in signals yi and y2 respectively. The MF-DFA plots Fq{s) 
vs s with fourth order detrending and q = —10, —8, —6, —4, —2, 2,4, 6,8, 10 
for signals yi and 2/2 are shown in Fig^l^c) and Fig^fd) respectively. For the 
monofractal data, the trends introduce spurious crossover at s ~ 2.2 in the log- 
log plot of the Fq{s) vs s for a given q. However, the nature of the log-log plots 
fail to show appreciable change with varying q, Fig^c). For multifractal data 
with trends, spurious crossovers are still noted at s ^ 2.2 in the log-log plot of 
the Fq{s) vs s for a given q. However, in this case, Fig^Jd), the nature of the 
log-log plots show a significant change with varying q indicating multifractal 
scaling in the given data, unlike the case with monofractal data with trends. 
Parameters: 

Ai = 6, A2 = 3, A3 = 2, Ti = 26, T2 = 24, T3 = 22,7Vi = 21^ B^ = 
6000, B2 = 3000, Tib = 2^, T2b = 2\ N2 = 2^^. 

The data sets are available from 
http://www.physionet.org/physiobank/database/synthetic/tns/ 
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Figure 4: MF-DFA studies of monofractal and multifractal data sets superim- 
posed with multiple sinusoidal trends. 
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